********************************************************************************************
**** QTREATTEST
**** This program compares means and equality of means at different quantiles
********************************************************************************************

	capt prog drop qtreattest

	program qtreattest, eclass

	syntax varlist [if] [in], [ * ]
	*marksample touse
	*markout `touse'
	tempname c_1 se_1 t_1 c_2 se_2 t_2 c_25 se_25 t_25 c_3 se_3 t_3 c_4 se_4 t_4 c_5 se_5 t_5 c_6 se_6 t_6 c_7 se_7 t_7 c_75 se_75 t_75 c_8 se_8 t_8 c_9 se_9 t_9 cf_1 sef_1 tf_1 cf_2 sef_2 tf_2 cf_25 sef_25 tf_25 cf_3 sef_3 tf_3 cf_4 sef_4 tf_4 cf_5 sef_5 tf_5 cf_6 sef_6 tf_6 cf_7 sef_7 tf_7 cf_75 sef_75 tf_75 cf_8 sef_8 tf_8 cf_9 sef_9 tf_9 obs
                           
	foreach var of local varlist {
	
	* BASELINE
			
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(10)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_1' = nullmat(`c_1'), b[1,1]
			mat `se_1' = nullmat(`se_1'), sqrt(v[1,1])
			mat `t_1' = nullmat(`t_1'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(20)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_2' = nullmat(`c_2'), b[1,1]
			mat `se_2' = nullmat(`se_2'), sqrt(v[1,1])
			mat `t_2' = nullmat(`t_2'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
		
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(25)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_25' = nullmat(`c_25'), b[1,1]
			mat `se_25' = nullmat(`se_25'), sqrt(v[1,1])
			mat `t_25' = nullmat(`t_25'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(30)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_3' = nullmat(`c_3'), b[1,1]
			mat `se_3' = nullmat(`se_3'), sqrt(v[1,1])
			mat `t_3' = nullmat(`t_3'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(40)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_4' = nullmat(`c_4'), b[1,1]
			mat `se_4' = nullmat(`se_4'), sqrt(v[1,1])
			mat `t_4' = nullmat(`t_4'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(50)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_5' = nullmat(`c_5'), b[1,1]
			mat `se_5' = nullmat(`se_5'), sqrt(v[1,1])
			mat `t_5' = nullmat(`t_5'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(60)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_6' = nullmat(`c_6'), b[1,1]
			mat `se_6' = nullmat(`se_6'), sqrt(v[1,1])
			mat `t_6' = nullmat(`t_6'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(70)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_7' = nullmat(`c_7'), b[1,1]
			mat `se_7' = nullmat(`se_7'), sqrt(v[1,1])
			mat `t_7' = nullmat(`t_7'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
		
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(75)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_75' = nullmat(`c_7'), b[1,1]
			mat `se_75' = nullmat(`se_7'), sqrt(v[1,1])
			mat `t_75' = nullmat(`t_7'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
		
		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(80)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_8' = nullmat(`c_8'), b[1,1]
			mat `se_8' = nullmat(`se_8'), sqrt(v[1,1])
			mat `t_8' = nullmat(`t_8'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)

		xi: qui qreg2 `var' mother $controls if round == 2010, cluster(id_mun) quantile(90)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `c_9' = nullmat(`c_9'), b[1,1]
			mat `se_9' = nullmat(`se_9'), sqrt(v[1,1])
			mat `t_9' = nullmat(`t_9'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)


	* FOLLOW-UP
		
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(10)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_1' = nullmat(`cf_1'), b[1,1]
			mat `sef_1' = nullmat(`sef_1'), sqrt(v[1,1])
			mat `tf_1' = nullmat(`tf_1'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(20)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_2' = nullmat(`cf_2'), b[1,1]
			mat `sef_2' = nullmat(`sef_2'), sqrt(v[1,1])
			mat `tf_2' = nullmat(`tf_2'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(25)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_25' = nullmat(`cf_25'), b[1,1]
			mat `sef_25' = nullmat(`sef_25'), sqrt(v[1,1])
			mat `tf_25' = nullmat(`tf_25'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(30)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_3' = nullmat(`cf_3'), b[1,1]
			mat `sef_3' = nullmat(`sef_3'), sqrt(v[1,1])
			mat `tf_3' = nullmat(`tf_3'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(40)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_4' = nullmat(`cf_4'), b[1,1]
			mat `sef_4' = nullmat(`sef_4'), sqrt(v[1,1])
			mat `tf_4' = nullmat(`tf_4'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(50)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_5' = nullmat(`cf_5'), b[1,1]
			mat `sef_5' = nullmat(`sef_5'), sqrt(v[1,1])
			mat `tf_5' = nullmat(`tf_5'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(60)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_6' = nullmat(`cf_6'), b[1,1]
			mat `sef_6' = nullmat(`sef_6'), sqrt(v[1,1])
			mat `tf_6' = nullmat(`tf_6'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
			 
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(70)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_7' = nullmat(`cf_7'), b[1,1]
			mat `sef_7' = nullmat(`sef_7'), sqrt(v[1,1])
			mat `tf_7' = nullmat(`tf_7'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
		
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(70)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_75' = nullmat(`cf_75'), b[1,1]
			mat `sef_75' = nullmat(`sef_75'), sqrt(v[1,1])
			mat `tf_75' = nullmat(`tf_75'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)
		
		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(80)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_8' = nullmat(`cf_8'), b[1,1]
			mat `sef_8' = nullmat(`sef_8'), sqrt(v[1,1])
			mat `tf_8' = nullmat(`tf_8'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)

		xi: qui qreg2 `var' mother $controls if round == 2012, cluster(id_mun) quantile(80)
		
			matrix b = e(b)
			matrix v = e(V)
			mat `cf_9' = nullmat(`cf_9'), b[1,1]
			mat `sef_9' = nullmat(`sef_9'), sqrt(v[1,1])
			mat `tf_9' = nullmat(`tf_9'), 2*ttail(e(df_r),abs(b[1,1]/sqrt(v[1,1])))
			mat `obs' = nullmat(`obs'), e(N)	
													
	}
		
	foreach mat in c_1 se_1 t_1 c_2 se_2 t_2 c_25 se_25 t_25 c_3 se_3 t_3 c_4 se_4 t_4 c_5 se_5 t_5 c_6 se_6 t_6 c_7 se_7 t_7 c_75 se_75 t_75 c_8 se_8 t_8 c_9 se_9 t_9 cf_1 sef_1 tf_1 cf_2 sef_2 tf_2 cf_25 sef_25 tf_25 cf_3 sef_3 tf_3 cf_4 sef_4 tf_4 cf_5 sef_5 tf_5 cf_6 sef_6 tf_6 cf_7 sef_7 tf_7 cf_75 sef_75 tf_75 cf_8 sef_8 tf_8 cf_9 sef_9 tf_9 {
		mat coln ``mat'' = `varlist'
	}

	tempname b V
	mat `b' = `c_1'*0
	mat `V' = `b''*`b'
	eret post `b' `V'
	eret local cmd "qtreattest"
	
	mat temp = `obs'
	eret scalar N1 =  temp[1,1] 
	eret scalar N2 =  temp[1,2]
	eret scalar N25 = temp[1,3]  
	eret scalar N3 =  temp[1,4] 
	eret scalar N4 =  temp[1,5] 
	eret scalar N5 =  temp[1,6] 
	eret scalar N6 =  temp[1,7] 
	eret scalar N7 =  temp[1,8]
	eret scalar N75 = temp[1,9]
	eret scalar N8 =  temp[1,10]
	eret scalar N9 =  temp[1,11] 
	
	eret scalar NF1 =  temp[1,12] 
	eret scalar NF2 =  temp[1,13]
	eret scalar NF25 = temp[1,14]  
	eret scalar NF3 =  temp[1,15] 
	eret scalar NF4 =  temp[1,16] 
	eret scalar NF5 =  temp[1,17] 
	eret scalar NF6 =  temp[1,18] 
	eret scalar NF7 =  temp[1,19]
	eret scalar NF75 = temp[1,20]
	eret scalar NF8 =  temp[1,21]
	eret scalar NF9 =  temp[1,22] 
	
	foreach mat in c_1 se_1 t_1 c_2 se_2 t_2 c_25 se_25 t_25 c_3 se_3 t_3 c_4 se_4 t_4 c_5 se_5 t_5 c_6 se_6 t_6 c_7 se_7 t_7 c_75 se_75 t_75 c_8 se_8 t_8 c_9 se_9 t_9 cf_1 sef_1 tf_1 cf_2 sef_2 tf_2 cf_25 sef_25 tf_25 cf_3 sef_3 tf_3 cf_4 sef_4 tf_4 cf_5 sef_5 tf_5 cf_6 sef_6 tf_6 cf_7 sef_7 tf_7 cf_75 sef_75 tf_75 cf_8 sef_8 tf_8 cf_9 sef_9 tf_9 obs {
		eret mat `mat' = ``mat''
	}
	
	end
